package com.greenrift.geneticshelper;

import android.app.Activity;
import android.os.Bundle;
import android.util.Log;
import android.view.View;
import android.view.Window;
import android.view.WindowManager;
import android.view.View.OnClickListener;
import android.widget.Button;
import android.widget.EditText;
import android.widget.TextView;
import android.widget.Toast;
import cern.jet.stat.Probability;

public class GeneticsHelperHWE extends Activity
{
  private int homo1 = 0;
  private int het = 0;
  private int homo2 = 0;
  private double result = -1d;

  private OnClickListener calculateClickListener = new OnClickListener()
  {
    public void onClick( View view )
    {
      EditText h1 = (EditText) findViewById( R.id.homozyg1_edit );
      EditText ht = (EditText) findViewById( R.id.het_edit );
      EditText h2 = (EditText) findViewById( R.id.homozyg2_edit );

      try
      {
        homo1 = Integer.parseInt( h1.getText().toString() );
      }
      catch( NumberFormatException ex )
      {
        Toast.makeText( getBaseContext(), "Homozygous Allele 1 needs to be an integer.", Toast.LENGTH_SHORT );
        return;
      }
      try
      {
        het = Integer.parseInt( ht.getText().toString() );
      }
      catch( NumberFormatException ex )
      {
        Toast.makeText( getBaseContext(), "Heterozygous needs to be an integer.", Toast.LENGTH_SHORT );
        return;
      }
      try
      {
        homo2 = Integer.parseInt( h2.getText().toString() );
      }
      catch( NumberFormatException ex )
      {
        Toast.makeText( getBaseContext(), "Homozygous Allele 2 needs to be an integer.", Toast.LENGTH_SHORT );
        return;
      }
      if( homo1 == 0 && het == 0 && homo2 == 0 )
      {
        Toast.makeText( getBaseContext(), "All values are 0.  No HWE can be calculated.", Toast.LENGTH_SHORT );
      }
      if( homo1 < 0 || het < 0 || homo2 < 0 )
      {
        Toast.makeText( getBaseContext(), "A value is < 0.  No HWE can be calculated.", Toast.LENGTH_SHORT );
      }
      int pop = homo1 + het + homo2;
      int apop = 2 * homo1 + het + het + 2 * homo2;
      float a1freq = (float) ( ( 2 * homo1 ) + het ) / (float) apop;
      float a2freq = 1 - a1freq;
      Log.d( "HWE", "FREQS = " + String.valueOf( a1freq ) + " : " + String.valueOf( a2freq ) );
      float p2 = a1freq * a1freq;
      float q2 = a2freq * a2freq;
      float pq = 2 * a1freq * a2freq;

      float expp2 = p2 * pop;
      float expq2 = q2 * pop;
      float exppq = pq * pop;

      TextView v = (TextView) findViewById( R.id.exp_rare_view );
      v.setText( String.valueOf( expp2 ) );
      v = (TextView) findViewById( R.id.exp_het_view );
      v.setText( String.valueOf( exppq ) );
      v = (TextView) findViewById( R.id.exp_common_view );
      v.setText( String.valueOf( expq2 ) );

      Log.d( "HWE", "EXP = " + String.valueOf( expp2 ) + " : " + String.valueOf( exppq ) + " : " + String.valueOf( expq2 ) );
      if( expq2 < 5 || expp2 < 5 || exppq < 5 )
      {
        // do exact
        Log.d( "HWE", "Do Exact" );
        result = exact( het, homo1, homo2 );

      }
      else
      {
        double chi1 = 0f;
        double chi2 = 0f;
        double chi3 = 0f;

        chi1 = ( ( (double) homo1 - expp2 ) * ( (double) homo1 - expp2 ) ) / expp2;
        chi2 = ( ( (double) homo2 - expq2 ) * ( (double) homo2 - expq2 ) ) / expq2;
        chi3 = ( ( (double) het - exppq ) * ( (double) het - exppq ) ) / exppq;
        Log.d( "HWE", "CHI = " + String.valueOf( chi1 ) + " : " + String.valueOf( chi2 ) + " : " + String.valueOf( chi3 ) );
        chi1 += chi2 + chi3;
        TextView vc = (TextView) findViewById( R.id.exp_chisq_view );
        vc.setText( String.valueOf( chi1 ) );

        // ChiSquare cs = new ChiSquare(1, RandomEngine.makeDefault());
        // result = cs.cdf(chi1);
        result = Probability.chiSquareComplemented( 1, chi1 );
      }

      if( result != -1d )
      {
        TextView pval = (TextView) findViewById( R.id.pval_view );
        pval.setText( Double.toString( result ) );
      }
    }
  };

  private double exact( int obs_hets, int obs_hom1, int obs_hom2 )
  {
    int obs_homc = obs_hom1 < obs_hom2 ? obs_hom2 : obs_hom1;
    int obs_homr = obs_hom1 < obs_hom2 ? obs_hom1 : obs_hom2;

    int rare_copies = 2 * obs_homr + obs_hets;
    int genotypes = obs_hets + obs_homc + obs_homr;
    if( genotypes <= 0 )
    {
      return -1;
    }

    double[] het_probs = new double[ rare_copies + 1 ];
    int i;
    for( i = 0; i <= rare_copies; i++ )
      het_probs[ i ] = 0.0;
    /* start at midpoint */
    int mid = rare_copies * ( 2 * genotypes - rare_copies ) / ( 2 * genotypes );
    /* check to ensure that midpoint and rare alleles have same parity */
    if( ( rare_copies & 1 ) == 1 ^ ( mid & 1 ) == 1 )
      mid++;

    int curr_hets = mid;
    int curr_homr = ( rare_copies - mid ) / 2;
    int curr_homc = genotypes - curr_hets - curr_homr;

    het_probs[ mid ] = 1.0;
    double sum = het_probs[ mid ];
    for( curr_hets = mid; curr_hets > 1; curr_hets -= 2 )
    {
      het_probs[ curr_hets - 2 ] = het_probs[ curr_hets ] * curr_hets * ( curr_hets - 1.0 )
          / ( 4.0 * ( curr_homr + 1.0 ) * ( curr_homc + 1.0 ) );
      sum += het_probs[ curr_hets - 2 ];

      /*
       * 2 fewer heterozygotes for next iteration -> add one rare, one common
       * homozygote
       */
      curr_homr++;
      curr_homc++;
    }

    curr_hets = mid;
    curr_homr = ( rare_copies - mid ) / 2;
    curr_homc = genotypes - curr_hets - curr_homr;
    for( curr_hets = mid; curr_hets <= rare_copies - 2; curr_hets += 2 )
    {
      het_probs[ curr_hets + 2 ] = het_probs[ curr_hets ] * 4.0 * curr_homr * curr_homc
          / ( ( curr_hets + 2.0 ) * ( curr_hets + 1.0 ) );
      sum += het_probs[ curr_hets + 2 ];
      /*
       * add 2 heterozygotes for next iteration -> subtract one rare, one common
       * homozygote
       */
      curr_homr--;
      curr_homc--;
    }

    for( i = 0; i <= rare_copies; i++ )
      het_probs[ i ] /= sum;

    /*
     * alternate p-value calculation for p_hi/p_lo double p_hi =
     * het_probs[obs_hets]; for (i = obs_hets + 1; i <= rare_copies; i++) p_hi
     * += het_probs[i];
     * 
     * double p_lo = het_probs[obs_hets]; for (i = obs_hets - 1; i >= 0; i--)
     * p_lo += het_probs[i];
     * 
     * double p_hi_lo = p_hi < p_lo ? 2.0 * p_hi : 2.0 * p_lo;
     */

    double p_hwe = 0.0;
    /* p-value calculation for p_hwe */
    for( i = 0; i <= rare_copies; i++ )
    {
      if( het_probs[ i ] > het_probs[ obs_hets ] )
        continue;
      p_hwe += het_probs[ i ];
    }
    p_hwe = p_hwe > 1.0 ? 1.0 : p_hwe;
    return p_hwe;

  }

  /** Called when the activity is first created. */
  @Override
  public void onStart()
  {
    super.onStart();
  }

  @Override
  public void onResume()
  {
    super.onResume();
  }

  @Override
  public void onPause()
  {
    super.onPause();
  }

  private void setFullscreen()
  {
    requestWindowFeature( Window.FEATURE_NO_TITLE );
    getWindow().setFlags( WindowManager.LayoutParams.FLAG_FULLSCREEN, WindowManager.LayoutParams.FLAG_FULLSCREEN );
  }

  @Override
  public void onCreate( Bundle savedInstanceState )
  {
    super.onCreate( savedInstanceState );

    setFullscreen();

    setContentView( R.layout.hwe );

    Button button = (Button) findViewById( R.id.calchwe_button );
    button.setOnClickListener( calculateClickListener );
  }
}
